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ABSTRACT 

A  simple  time  series  model  for  bivariate  exponential  variables 
having  first-order  autoregressive  structure  is  presented.  The  linear 
random  coefficient  difference  equation  model  is  an  adaptation  of  the  New 
Exponential  Autoregressive  model  (NEAR  (2)).  The  process  is  Markovian 
in  the  bivariate  sense  and  has  correlation  structure  analogous  to  that 
of  the  Gaussian  AR(1)  bivariate  time  series  model.  The  model  exhibits  a 
full  range  of  positive  correlations  and  cross-correlations.  With  some 
modification  in  either  the  Innovation  or  the  random  coefficients,  the 
model  admits  some  negative  values  for  the  cross-correlations.  The 
marginal  processes  are  shown  to  have  correlation  structure  of  ARMA(2,1) 
models. 


Key  Words:  Time  Series;  Bivariate  Exponential  Distribution; 

Autoregressive  Models;  NEAR(2);  ARMA(2,1)  Models;  Gaussian 
AR(1)  Bivariate  Time  Series  Model. 


1.  INTRODUCTION 


The  homogeneous  Poisson  process  is  a  basic  model  for  point  processes 
(series  of  events)  and  can  be  characterized  as  a  process  in  which  the 
intervals  between  events  are  independent  and  identically  exponentially 
distributed  random  variables.  One  generalization  of  the  homogeneous 
Poisson  process  has  been  to  relax  the  condition  on  independence  of 
intervals.  The  earliest  attempts  to  do  this  go  back  to  Wold  (1948)  and 
Cox  (1955).  See  Cox  and  Lewis  (1966,  Ch.  7)  for  details.  Later,  this 


11/5/86 


LEWIS/DEWALD 


y  Codes 

)..d/or 

cial 


rsan 


need  for  point  processes  with  correlated  but  marginally  exponential 
distributions  inspired  the  creation  of  a  class  of  ARMA  time  series 
models  for  the  univariate  exponential  processes  analogous  to  that  for 
the  Gaussian  processes  (Caver  and  Lewis,  1980;  Lawrance  and  Lewis,  1981, 
1985).  Beyond  their  utility  in  point  processes,  these  models  are 
utilized  in  studying  queues  (Jacobs,  1978;  1980),  inventory  and  water 
resources  problems,  and  other  situations  where  non-negative  random 
variables  are  appropriate  for  Inputs.  Using  antithetic  variables  as 
innovation  or  error  structure,  these  models  also  allow  for  negative 
dependence  (see  Lewis,  1986,  for  a  survey). 

We  demonstrate  that  one  of  the  more  recent  exponential  models  with 
ARMA  correlation  strucutre,  NEAR(2),  by  Lawrance  and  Lewis  (1985),  lends 
itself  naturally  to  a  first-order  bivariate  autoregressive  process  which 
is  Markovian  in  the  bivariate  sense,  and  has  a  structure  analogous  to 
that  of  the  Gaussian  AR(1)  bivariate  time  series  model  as  given  in  Tiao 
and  Box  (1981),  and  elsewhere.  Early  attempts  to  construct  bivariate 
exponential  time  series  and  bivariate  Poisson  processes  are  discussed  in 
Cox  and  Lewis  (1971).  The  present  model  is  broader  and  simpler  than  any 
previously  obtained,  including  the  MEAR(l)  model  by  Raftery  (1982), 
which  turns  out  to  be  a  special  case  of  the  present  model  for  the 
bivariate  case,  and  the  models  of  Jacobs  (1978,  1980),  Lewis  and  Shedler 
(1977),  and  Gaver  and  Lewis  (1980).  All  of  these  models  were  defined 
using  the  NEAR(l)  structure,  and  this  does  not  allow  for  the  breadth 
obtained  by  using  the  NEAR(2)  structure. 

With  the  general  model  it  is,  for  example,  possible  to  construct 
simple  bivariate,  serially  correlated  models  for  the  successive  up  and 
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down  times  which  occur  in  reliability  studies,  or  for  correlated  and 
cross-correlated  service  and  Interarrival  sequences  In  queueing  systems. 

2.  THE  UNIVARIATE  NEAR (2}  MODEL 

By  the  NEAR(2)  theorem  (Lawrance  and  Lewis,  1985),  we  know  that  two 
(possibly  dependent)  random  variables  with  marginally  Exponential (X) 
distributions  can  be  combined  with  three  scaled  (possibly  dependent) 
exponentials  which  are  Independent  of  the  first  pair,  to  give  another 
random  variable  that  is  marginally  exponentially  distributed  with 
parameter  X.  We  assume  in  the  remainder  of  this  presentation  that 
X  =  1.  The  essence  of  the  NEAR(2)  theorem  is  repeated  without  proof. 

THEOREM  (Lawrance  and  Lewis) 

If  (En)  is  an  i.i.d.  sequence  of  exponential  random  variables  and  the 
parameters  a^,  0^,  a Pg  are  suc^  t*iat  ®  <  aj  <1  for  1=1,2; 


<  al  +  a2  <  1 

and  0  <  0J,  P2 

<  1 

and  if 

IS 

x° 

Eq  ,  then 

'  P.Xn-l 

w.p.  <Xj 

E 

n 

w.p. 

1-P2'p3 

11 

*C 

P2Xn-2 

w.p.  a2 

+ 

b2En 

w.p. 

p2 

0 

w.p.  l-ctj-ag 

,b3E„ 

w.p. 

P3, 

defines  a  stationary  sequence  (Xn)  which  is  marginally  exponentially 
distributed  with  X=l.  The  quantities  bg,  b^,  Pg.  P3  have  values  in 
the  interval  (0,1)  and  are  defined  in  Lawrance  and  Lewis  (1985);  the 
mixture  of  the  exponential  random  variables  on  the  right  in  (2.1)  is 
called  the  innovation  random  variable. 
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A  ■pedal  case  of  the  NEAR(2)  aodel  which  la  not  covered  by  the 
theorem,  but  obtained  directly  from  the  usual  analysis  with  transforms, 
is  the  case  when  0^  =  0,  =  1.  In  this  case,  still  assuming 

°1  +a2  *  we  ^iave 


X  = 
n 


n-1 


n-2 


w.p. 

w.p.  a2 

W.p.  l-Oj-Og 


+  (l-o1-a2)En.  (2.2) 


This  is  the  TEAR(2)  model.  It  has  two  advantages  for  our  purposes;  one 
that  It  only  uses  two  parameters  instead  of  four  in  the  complete  NEAR(2) 
model;  and  the  other  that  it  has  a  simpler  innovation  structure.  The 
model  exhibits,  however,  a  "runs-up"  behavior  in  the  sample  paths. 
Nevertheless,  we  will  use  the  TEAR(2)  model  as  the  base  for  construction 
of  the  bivariate  exponential  model,  even  though  other  three-  and 
two-parameter  subclasses  of  the  NEAR(2)  process  may  be 
phenomenologically  more  suitable.  For  example,  in  hydrology  the  time 
series  which  occur  generally  exhibit  a  ”runs-down"  behavior.  However, 
the  ideas  in  developing  the  bivariate  process  are  the  same. 


3.  BIVARIATE.  MARKOVIAN  TIME  SERIES  FOR  EXPONENTIAL  RANDOM  VARIABLES 
Let  {^n»Yn)  be  a  bivariate  sequence  of  random  variables,  such  that 
(X^)  and  {Y^}  are  each  stationary  sequences  of  marginally 
exponentially  distributed  random  variables  with  X=l.  Then  we  define 
our  first  bivariate  exponential  model,  using  (2.2),  as 

X  =  k[^X  ,  +  kS"5Y  .  +  (l-a11-ain)E  (3.1) 

n  11  n-1  12  n-1  v  11  12'  n  '  ' 
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and 


\  -  4?’  Xn-1  *  <3'2> 

In  this  equation 

(0.1)  w.p.  aJ2 

(1.0)  w.p.  an 

(0.0)  w.p.  l-Ojj-otjg 

and  are  serially  Independent; 

(0.1)  w.p. 

K^}  =  (1.0)  w.p.  02J 

(0.0)  w.p.  x_a2i_a22 

are  also  serially  independent  and  Independent  of  (K^.  4*^)  *or 
n.  Likewise  (En)  and  {E^}  are  Independent  Innovation  sequences  of 
l.l.d.  Exponential (A=l)  random  variables.  Finally,  we  Insist  that 

aH  +  o12  <  1  and  «22  +  “2!  <  *•  If  «n  +  “i2  =  1 

°22  +  “21  =  * ’  t*ien  t^le  process  Is  not  ergodic.  In  this  case, 

is  always  one  of  the  pairs  (Ej.Ej).  (Ej.Ej),  (Ej.Ej)  (Ej.Ej)  for  all 

n.  Note  that  ^n_j  and  ^n_l  in>  ®-8*»  Equation  (3.1)  are  not 

independent;  this  is  the  reason  for  the  use  of  the  NEAR(2)  construction. 

We  observe  that  (Xn) ,  (Yn)  are  each  TEAR(2)  constructions  in  that 

Yn_j  in  (3.1)  replaced  ^n_j  in  (2.1),  and  similarly  Xn_j  in  (3.2) 

replaced  Y  n.  Also,  (X  ,Y  )  is,  given  the  value  of  (X  .  ,Y  , ), 

n-z  n  n  n-i  n-i 

completely  Independent  of  previous  values  (X^.Y^)  *or  ^  £  n-2 
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all  n.  Again  *(Xn|Xn_j  =  x-  Ynl  =  y)  is  a  linear  function  of  (x,y) 

and  likewise  for  Y  .  Thus,  the  model  is  Markovian  in  the  bivariate 

n 

sense,  both  structurally  and  in  the  sense  of  expectations. 

Note,  however,  that  {Xr}  and  {Yr}  are  Qfit  marginally  TEAR(l) 
processes,  and,  in  fact,  are  not  marginally  Markovian.  In  particular, 
their  correlations  are  not  necessarily  geometrically  decaying,  as  will 
be  demonstrated  later.  If  =  a99  =  0,  we  obtain  the  bivariate 

exponential  process  described  by  Gaver  and  Lewis  (1980).  The  process 
(Yn)  is  then  called  the  dual  of  the  process  (Xn). 

Clearly,  other  special  combinations  of  the  give  other  models 

with  1.  2.  or  3  parameters.  Of  course,  there  could  be  from  5  to  8  total 
parameters  by  using  the  more  general  NEAR(2)  construction.  Thus,  from 
(2.1),  (3.1),  and  (3.2),  replace  K by  0^  K^.  0  <  0^  <  1.  One 
requires  then,  however,  the  more  complicated  innovation  structure  given 
in  (2.1)  for  the  innovation  in  (3.1)  and  (3.2). 

4.  CORRELATION  STRUCTURE  OF  THE  MODEL 

The  marginal  processes  defined  above  are  exponentially  distributed 

2 

with  mean  1/A  and  variance  1/A  .  Thus,  if  we  set  the  scale  parameter 
A  equal  to  1,  correlations  and  covariance  are  equal.  Since  we  assume 
(Yn).  (Xn)  are  Jointly  stationary,  we  have 

7^(0  =  *(XnXn_p  -  1.  I  =  0.  1,  2.  •••.  (4.1) 

=  *(YnYn_p  -  1.  i  =  0.  1.  2.  •••.  (4.2) 

for  the  autocovariance  functions.  Similarly, 
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(4.3) 


V<‘>  *  *<X„V«>  -  1  *  =  <>•*>•*  *. 

is  the  cross-covariance  function. 

The  following  results  are  immediately  obtained  using  (3.1).  (3.2), 
(4.1)  -  (4.3)  and  the  stationary  assumption: 

^(0)  =  *rYY(0)  =  1:  (44) 

nrXY(0)  =  =  {a22a12  +  all“21)/{1  ~  (a22all  +  “21a12)}:  (4  5) 


^.,.(0  =  f°r 

(4.6) 

'ryy(i)  =  TyvC-^)  for  all  t ; 

(4.7) 

Txy(^)  =  *or 

(4.8a) 

'Tyx(^)  =  for  all 

(4.8b) 

Autocovariance  functons  are  even  functions  of  the  lag  i,  but  the 
cross-covariance  is  not. 

The  range  of  ^^(0)  =  as  S*ven  *n  (4.5)  is  the  full  range 
of  non-negative  values  up  to  one.  Using  the  Generalized  Reduced 
Gradient  (GRG)  algorithm,  we  find  that  the  maximum  takes  place  along  the 
boundary  ctjj  +  a ^  +  e  =  1  and  +  “21  +  £  =  1  for  arbitrary  e  >  0, 
in  which  case  we  obtained  covariances  as  close  to  one  as  desired.  Of 
course,  for  strictly  positive  a^j  that  satisfy  the  original 
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constraints,  **^(0)  cannot  be  negative.  Ve  show  In  Section  5  how  to 
Incorporate  negative  correlations  into  the  model. 

Using  Equations  (3.1),  (3.2),  (4.1)-(4.8)  and  the  stationary  assumption, 
we  obtain  the  following  recursion  equations  for  ^XX^’  7YY^‘ 

*  7YX^  : 


=  “iTxx^-^ 

+  a127YX^-1^’ 

t  =  1.2.  —  . 

(4.9) 

=  a22'TYY^-1^ 

4  a2l'rXY^~1^* 

t  *  1,2, ••• , 

(4.10) 

4  “227YX(I_1>* 

i  =  1,2,  —  , 

(4.11) 

=  “ig^YY^-15 

+  an'rxY^~^^  * 

l  =  1,2,  —  . 

(4.12) 

It  remains  to  show  that  these  equations  correspond  to  those  given 
by,  for  example,  Tiao  and  Box  (1981),  for  the  AR(1)  Gaussian  bivariate 
time  series.  We  introduce  the  following  matrix  notation: 


r(i) 


7XX^  yXY^ 

7yx^^ 


i  =  0.1.2.  —  .  (4.13) 


is  the  covariance  matrix  for  lag  i.  Equations  (3.1)  and  (3.2)  become 


X  =  K  X  . 
n  n  n-1 


+  AE 


n 


(4.14) 
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where 


1_“ll“®12  0 

Xn 

= 

:  X  * 

0  1_a22_a21 

n 

Y 

n 

and 


K  = 
n 


E  = 
n 


n 


E* 

n 


K(n) 

K(n> 

K11 

*12 

Kin> 

*22 

The  matrix-valued  random  variable  K  assumes  9  different 

n 

matrix-values  independent  of  and  En  and  retains  the  property  that 

a  for  all  1,  J.  n.  Then 


<(Xn)  =  «(K„)«(Xn_i)  +  A  <(En)  (4.15) 


“ll  a12 

'  1  ' 

'  1_an"°i2  0 

1  ' 

'  1 

= 

+ 

= 

“21  “22 

1 

0  1'a22'°21 

1 

»  « 

1 

Also. 

r(°)  =  «(x$  -  «(xn)«(xj)  =  «(xnxj)  -  [!  | 

■  ‘<Wn-l  *  “»>T>  -  [1  !] 

■  +  ‘^(“n)7)  '  [!  !] 
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■  *«KnXn-l  ♦  “„>(Vn-l>T>  +  ««K„X„-1  *  “»><")/>  “  [l  i] 
=  «<K„Xn-l  (tt-ll’l  +  ,(*E„(Kr>Vl)T>  +  ‘<KnXn-l  (^>T> 

+  n*EX*T)-[!  }]• 

Thus 

no)  .  «(K„xn-ixli  <>  ♦  *  >(^y> 

*  ,(tEnEy>  -[!  1]  •  (•‘■16) 

Solving  (4.16)  yields  the  results  given  by  Equations  (4.4)  and  (4.5). 
Furthermore 

no  =  -  «oy«oi>  =  ‘oyO  -  [1  1] 

-  *«KA-1  +  Wll*  -  [!  1] 

=  ,(Kn>r<°>  +  *<K„>  [1  1]  +  A«EnXll>  ‘  [l  }]• 

Thus 


r(i)  =  i(Kn)r(o)  = 


°n 

°21 


“12 

°22 


no). 


(4. 17) 
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Again 


r(-i)  = 

'<xX.i>  - 

[1  1] 

= 

‘<X„<Xn*lXn 

*  -  [1  !] 

= 

*<x„xX+i> 

♦‘<x„eLi>aT-  [1  !]• 

Thus 

r(- 

-1)  =  r(0)i(Kj+1)  =  r(i)T 

(4.18) 

Clearly  T(- 

-i)T  =  rci). 

Because  of  (4.8a,b)  and  (4.13), 

we  have 

n-t)  ■=  nty 

(4.19) 

and 


r(0 


an 

“21 


a12 

°22 


e 

r(0) 


for  all  t  =  1,2, •••. 


(4.20) 


Nov  Eolation  (4.20)  is  the  result  given  by  Tiao  and  Box,  (1981, 
p  8C *).  Therefore,  the  auto-  and  cross-covariances  in  this  bivariate 
exponential  process  in  general  decay  gradually  to  zero  as  |l| 
Increases,  for  the  values  of  that  satisfy  the  conditions  in 

Section  3.  Numerical  examples  are  given  in  Section  6,  and  further 
discussion  of  the  auto-correlations  is  presented  in  Section  7. 
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In  principle,  one  can  use  standard  estimates  of  auto-  and 
cross-correlations  (Kendall  and  Stuart.  1963,  pp.  661-662)  to  obtain 
estimates  of  T(0)  and  r(l)  as 


and 


C(0) 


1 

W0* 


i 


(4.21) 


C(l) 


(^(1)  <^0) 

Cyx(l) 


(4.22) 


Then  (4.20)  yields  estimates  of  a^,  a^j.  a^,  as 


all 

a12 

“21 

A 

°22 

(C(l)HC(O)) 


(4.23) 


Like  most  moment  estimates,  these  are  likely  to  be  fairly  poor  for  small 
sample  size. 

An  example  Is  given  in  Table  4.1  where  we  have  performed  a 
simulation  of  the  estimation  of  a^,  a^.  and  a99,  based  on  (4.23). 
Here  10  replications  were  used  and  the  true  values  were  =  0.70; 
otj2  =  0.20;  Ogj  =  0.30;  =  0.50.  The  10  series  were  run  out  to 
n  =  600  and  the  first  100  values  were  discarded  as  being  nonstationary. 
The  boxplots  In  the  tables  show  the  distribution  of  the  tabulated  values 
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In  the  columns .  The  estmates  show  no  evidence  of  departure  from 

normality,  although  the  number  of  replications  Is  small. 


TABLE  4.1 


Replication  =  0.70  =  0.20  =  0.30  ap?  =  0.50 


JL 

s 

* 


.«J 


.70 


.17 


.71 


.71 


.10 


.74 


.11 


.*• 


0 


.si 


.17 


.to 


.17 


.11 


.11 


.11 


.17 


.11 


.10 


.10 


.47 


.11 


.41 


.10 


.10 


o 

* 

T 


.10 


.11 


.11 


.11 


.17 


.41 


.10 


.10 


.00 


10 


.70 


.10 


.11 


.10 


5.  NEGATIVE  DEPENDENCE 

As  pointed  out  In  the  last  section,  the  first  bivariate  exponential 

model  discussed  above  can  produce  only  positive  correlations.  Negative 

correlation  is  obtainable  in  one  of  two  ways  -  through  correlation  in 

the  bivariate  Innovation  {E^.F^}  and/or  through  the  correlation  in  the 

attenuation  controlled  by  the  K^.  These  ideas  are  relaxations  of  the 

conditions  established  in  Section  3  where  {E  }  and  {E*}  are  assumed 

1  n'  1  n' 


11/5/86 


13 


LEWIS/DEWALD 


to  be  i.i.d.  exponential  sequences  and  (K^,  K^}  18  811  independent 

identically  distributed  bivariate  sequence  Independent  of  {K^, 
for  all  n. 

It  is  easy  to  see  that  these  schemes  affect  directly  only  the  zero 
lag  covariances,  r(0),  and  that  (4.20)  still  holds.  Thus,  it  is  only 
necessary  here  to  derive  the  effect  on  T(0)  of  correlation  in  the 
innovation  or  attenuation  sequences. 

5.i.  Correlated  Innovations. 

As  a  first  example,  consider  the  model  of  Section  3  and  let 


°"<VW 


a  k=0 
0  otherwise. 


(5.1) 


Now  Equation  (4.5)  becomes 


’XY<°> 


^I1^21+a12a22^  +  ^1-an-0i2^1-a2l"a22^CT 
1  '  (“ll“22Wlftl>  ' 


Moran  (1967)  showed  that  if  (En)  and  (E^)  are  pairs  of  exponential 

random  variables  with  A  =  1,  then  -.6449  £  a  $  1.  If  E'  =  a(E  ),  the 

n  n 

antithetic  transform  of  Er  (see  for  example  Gaver  and  Lewis,  1960), 
then  a  takes  on  its  maximum  negative  value.  It  is  clear  that  if 
a  <  0,  then  (5.2)  could  be  negative  for  some  choices  of  the  a.  *s. 
The  process  cannot  be  uncoupled  as  described  previously. 

Estimation  of  a  proceeds  as  before.  The  parameters  a^,  a^,., 
,j.  a*. 2  are  estimated  by  (4.23)  and  then  (5.2)  can  be  easily  solved 
for  an  estimate  of  a. 
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5.11.  Correlation  In  the  Rand  oca  Coefficients. 

Consider  now  the  model  of  Section  3  In  the  form  of  (4.14).  Kr 
defines  a  quadruple  where  each  has  a  marginal  distribution  and  each 
row  has  an  Independent  bivariate  distribution  for  each  n.  If  Kr  is 
now  defined  so  that  the  rows  are  not  stochastically  Independent  of  each 
other,  then  (4.5)  becomes 


,(0)  = 


{  (al  la21*f<>12iQe22^  +  20 1 1.21+^12 . 22+trl  1 , 22+a12 . 21 


"  all,22  "  a12,21 


(5.3) 


where 


aiJ .in 


Cov(K[^,Kjj|+k^)  k  =  0 

0  otherwise. 


(5.4) 


This  correlation  must  be  introduced  In  such  a  way  that  the  marginal 
distributions  are  not  altered.  This  may  not  be  easy. 

The  following  example  using  for  all  n  demonstrates 

the  above  ideas.  Define  I(_  as  follows: 

W,P#  aH’ 
w.p.  a^ajg. 

w.p.  a12(1”a2i)'  (5.5) 

w.p.  1  -  (a!  1+a12’Kl21_a12a21  ^ * 
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where  a^,  a^.  <*21  “  h®**01*®*  and  also  0  £  all'*:,12+a21~a12a21  *  ** 

It  is  easily  verified  that  the  Marginal  distributions  of  the  quadruple 

are  retained  and,  specifically,  that 


«v  = 


“11 

.“21 


For  Kr  defined  by  (5.5)  we  have,  from  (5.3) 


''XY(0) 


gll(1~;ttllw12*a2l)> 

1  -  (“ll+a12<121) 


(5.6) 


(5.7) 


ThU  folio™  boom...  an  2l  .  -a,,*,,:  »12-22  =  -«,2«ni  a„  j 2  ■= 

aH (!““!!)  and  ai2,2\  =  °*  with  0  *  “ll4®^  *  1:  0  *  *  1; 

0  £  an+ai24a2l“°l2a21  ^  * '  Using  again  GRG,  we  find  ^^(0)  attains 
the  full  range  of  positive  correlations  as  before.  It  has  a  max 
negative  correlation  of  -.125  at  =  .25,  “  *5.  This 

corresponds  to  one  of  the  cases  where  /  0  for  any  n. 

Considering  the  bivariate  model  with  given  by  (5.5)  and  with 

all  =  we  see  t*iat  =  0  =  7^(2k) ,  k  =  0,1,2,»»».  Resolving 

(4.9)  -  (4.12)  yields 


and 


F(k)  = 


ik 


12 


“21 


k  =  1.2, 


P(0) 


0 

1 
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This  is  the  correlation  structure  of  the  dual  process  of  Gaver  and  Lewis 
(1980).  No  negative  correlation  is  possible  in  this  configuration. 

Lewis  and  Lawrance  (1961)  considered  the  model  with  correlation  in 
the  random  coefficients  and  made  and  kJ^'  have  the  maximum 

possible  negative  correlation,  obtaining 

^<°>  - 


where 


a  =  Cov(kM.  K<">)  .  c12  2l  <  0. 

Considering  the  bivariate  model  with  given  by  (5.5)  with 

m  =  0.  then  ^^(0)  *  and  the  resolution  of  (4.9)  -  (4.12) 

with  (5.7)  yields  the  following  structure.  Marginally  (Xn)  and  (Yn) 

are  univariate  TEAR(l)  processes,  since  the  Yn’s  do  not  appear  in  the 

definition  of  the  X  *s  and  vice  versa.  Moreover,  ttvv(#)  =  aft* 

n  at  1 1 

=  ^yx^’  e  Therefore,  there  is  no  negative  correlation 

nor  any  way  to  uncouple  the  bivariate  process  into  two  Independent 
univariate  TEAR(l)  processes  unless  =  0.  Of  course,  if  =  0, 
there  is  then  no  dependence  in  the  marginal  processes  either.  The 
bivariate  process  with  jt  0  is  shown  in  the  next  section  to  have 

Moran's  Bivariate  Exponential  distribution. 

Finally,  we  observe  that  if  an+ai2+a21  =  *  in  (5.7),  we  have 
*^(0)  =  0.  Thus,  the  pairs  are  uncorrelated  for  all  n  and 

yet  the  processes  (Xn)  and  (Yn)  are  not  at  all  Independent.  This 
can  be  seen  from  the  fact  that  -1^(1)  *  ctj g  /  0  and  *1^(1)  =  “21  ** 
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5.111.  Correlated  Innovations  and  Randoa  Coefficients. 

If  there  Is  dependence  from  both  the  Innovation  and  cross-coupled 
attenuation,  then  a  general  expression  for  (4.5)  becomes 


XY 


(0) 


{ (a2l“l l+flc22a12^^2crll , 21+2p12, 22^11 ,22^12. 21^+^1~0tll‘<r12^1“a21~a22^<y^ 


1  "  ^a12a21+a22all^  “  12.21^11. 22} 


(5.8) 


The  Interesting  point  to  be  made  here  Is  that  the  effects  are  additive 
In  the  numerator  of  (5.8).  This  general  model  is  called  the  HEAR(l) 
model. 

6.  NUMERICAL  EXAMPLES  OF  COVARIANCE  MATRICES. 

The  form  of  T(t)  given  in  (4.20)  can  be  misleading.  Since  the 
form  of  r(£)  Is  analogous  to  that  of  a  univariate  AR(1)  process,  one 
might  think  that  the  elements  of  T(A)  should  experience  a  geometric 
decay  as  |<|  Increases.  In  fact,  auto-  and  cross-covariances  can 
exhibit  many  different  patterns  depending  on  the  innovation  and  the 
attenuation  structures  of  the  model.  Some  Illustrative  examples  follow. 

In  Table  6.1  an  example  Is  described  where  ^xY^  I8  given  by 
(4.5).  Although  the  auto-covariance  functions  decrease  monotonlcally, 
it  Is  not  at  all  like  a  geometric  decay. 

In  Table  6.2,  an  example  corresponding  to  ^XY^^  as  given  in  (5.2) 
is  exhibited.  This  case  demonstrates  that  negative  correlation  in  the 
bivariate  Innovation  can  create  negative  cross-covariances  at  lag  zero. 
It  can  also  cause  oscillation  in  the  values  of  the  auto-  and  cross- 
covariance  functions  at  subsequent  lags.  Also,  note  for  this  example 
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that  simply  by  Interchanging  for  and  «21  *or  a22‘  w®  0811 
create  another  bivariate  series  with  the  same  T( 0),  but  with  auto¬ 
covariance  functions  '^(l)  and  t^at  are  a*most  Eaometrically 
decaying. 

In  Table  6.3.  we  give  a  final  example  corresponding  to  from 
(5.3).  It  is,  in  fact,  the  example  given  for  in  (5.5).  Again, 
tj£y(0)  is  negative.  However,  note  that  although  the  auto-covariance 
functions  are  generally  decreasing ,  the  movement  is  not  monotone  as  in 
Table  6.1,  and  certainly  not  decreasing  geometrically. 

It  is  apparent  from  all  examples  that  as  |f|  Increases.  T(<) 
approaches  the  zero  matrix.  In  fact,  since  the  eigenvalues  of  the 
matrix  are  less  than  one,  the  components  of  f(f)  all  decrease  to  zero. 


TABLE  6.1  Covariance  Matrix  Corresponding  to  (4.5) 


t 

“n  *  2 

°12  *  1 

-'„(«) 

“21  *  3 

’»<*> 

°22  =  5 

0 

1.0000 

0.5942 

0.5942 

1.0000 

1 

0.6159 

0.5971 

0.8188 

0.6783 

2 

0.5412 

0.4833 

0.6386 

0.5848 

3 

0.4466 

0.4040 

0.5371 

0.4840 

4 

0.3721 

0.3360 

0.4462 

0.4031 

5 

0.3096 

0.2796 

0.3714 

0.3354 

6 

0.2577 

0.2327 

0.3091 

0.2791 

7 

0.2144 

0.1936 

0.2572 

0.2323 

8 

0.1784 

0.1611 

0.2140 

0.1933 

9 

0.1485 

0.1341 

0.1781 

0.1609 

10 

0.1236 

0.1116 

0.1482 

0.1339 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

• 

25 

0.0089 

0.0081 

0.0106 

0.0097 
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TABLE  6.2  Covariance  Matrix  Corresponding  to  (5.2) 

otj  j  ■  .05 

a12  =  60  ®21  *  ,4° 

=  .01 

a  -  -.6449 

*  ’xx^ 

W**  7XY^ 

0 

1.0000 

-0.1411 

-0.1411 

1.0000 

1 

-0.0347 

0.3986 

0.5929 

-0.0464 

2 

0.2374 

-0.0099 

0.0018 

0.2367 

3 

0.0059 

0.0949 

0.1421 

0.0031 

4 

0.0572 

0.0033 

0.0090 

0.0569 

5 

0.0048 

0.0229 

0.0346 

0.0042 

6 

0.0140 

0.0021 

0.0043 

0.0139 

7 

0.0020 

0.0056 

0.0086 

0.0019 

8 

0.0035 

0.0009 

0.0016 

0.0035 

9 

0.0007 

0.0014 

0.0022 

0.0007 

10 

0.0009 

0.0002 

0.0005 

0.0009 

TABLE  6.3  Covariance  Matrix  Corresponding  to 

(5.3) 

ctj  i  s  .25  a\2  ~  = 

“22  “  25 

CT11.22  =  "-1  a12.22  =  "15  *11.22  “* 1875 

a12.2l“  0 

1  7xx(<)  'nW  'rxY(<) 

^yy^) 

0 

1.0000 

-0.1225 

-0.1225 

1.0000 

1 

0.1765 

0.3694 

0.5694 

0.2010 

2 

0.2659 

0.1629 

0.2629 

0.2780 

3 

0.1642 

0.1470 

0.2325 

0.1747 

4 

0.1293 

0.1024 

0.1629 

0.1367 

5 

0.0938 

0.0773 

0.1227 

0.0993 

6 

0.0696 

0.0568 

0.0903 

0.0739 

7 

0.0516 

0.0421 

0.0669 

0.0546 

8 

0.0382 

0.0312 

0.0495 

0.0404 

9 

0.0282 

0.0231 

0.0366 

0.0299 

10 

0.0209 

0.0171 

0.0271 

0.0221 
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7.  THE  NARGIKAL  PROCESS. 


Writing  K  -  l(iy 
relationship 


“ll 

“21 


a12 

“22 


we  have  from  (4.20)  the 


r(«)  -  K*r(0). 


(7.1) 


Further,  since  K  Is  a  square  matrix,  the  well-known  Qaylay- 
Hamilton  Theorem  proivdes  that  there  exists  constants  and 

such  that 


K2  -  K  -  *2I  =  0. 

Applying  this  result  to  (7.1)  yields 

r(i )  -  -  42r(«-2)  =0.  tl  2.  (7.2) 

The  difference  Equation  (7.2)  applies  to  each  element  of  the  matrix 

r(£).  Hence,  for  example,  the  autocorrelation  function  of  (X  } 

n 

satisfies 


PXX^  "  ♦lpkX^“1^  ”  *2PXX^“2)  s  1  2. 

Thus,  In  general,  (XR)  has  the  correlation  structure  of  an  ARMA(2,1) 
process.  In  addition,  the  lag  1  autocorrelation  may  be  derived  from 
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(7.1)  when  4*1.  In  particular.  'fjQfO)  “  “ll^XX^^  +  a12nrYX^^’ 
since  ^^(0)  =  ''fyyC®)  =  1.  w®  find  that  Pjq^O)  =  an  +  al2pYX^‘ 

Further,  the  constants  ^  and  ^  may  be  obtained  in  terns  of  the 
elements  of  K.  It  may  be  shown  directly  that  ^  =  (ctj.  op)  and 
-t^ct^).  Thus,  in  summary,  the  autocorrelation  functions 
for  the  marginal  processes  {XR}.  {Yn)  satisfy 

p(4)  =  (ttjj+a^Jpff-l)  +  (a12a21-ana22)p(e-2) .  (7.3) 

for  £2  2,  and 

PXX^  =  all  *  a12pYX^*  C7'4) 

Pyyt1)  =  a22  *  a21pXY^'  (7-5) 

In  general,  {X^}  and  (Yn)  have  correlation  structures  of 
ARMA(2.1)  processes.  In  any  particular  case,  their  actual  forn  will 
depend  upon  the  values  {a^ ^ } .  The  marginal  processes  will  liave 
ARMA(p.q)  correlations,  where  p  £  2  and  q  £  1.  Ve  illustrate  these 
ideas  with  some  examples. 

(1)  Using  a^2  =  =  0  results  in  the  X  and  Y  processes  being 

independent  TEAR(l)  processes.  From  (7.4)  and  (7.5), 

Pjq^(I)  =  “ir  PYY^^  =  a22’  an<^  usinS  (7.3)  w®  obtain 
pXx(*)  =  “jj  and  Pyy(4)  =  a^.  The  bivariate  system  has 
become 

X  =  k[^X  ,  +  (l-a.JE 
n  11  n-1  '  11'  n. 
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Y 

n 


K&l’Vl  *  “-“22>E;- 


where  P(k{“J-1)  =  «n  and  P(K^=1)  =  o^.  and  and 

K^)  are  Independent,  as  are  Er  and  E\ 

(11)  Using  =  0^2  =  0  leads  to  coupled  processes,  each  with 

AR(1)  correlation  structure.  Now,  =  a\2  and 

Pyytl)  =  “2jP*  where  p  =  p^(0).  From  (7.3)  we  may  derive 


PxX(2k)  =  (a12“21)k 


PxX(2k+1)  =  (a12“21)ka12P  J 


k  =  0, !,•••. 


This  Is  the  formulation  used  by  Caver  and  Lewis  (1980)  to 
derive  negatively  correlated  processes.  In  effect,  they  chose 

°12  =  “21  =  “*  ar'd  80 

k  , 
a  k  even 

a^p  k  odd. 

Now,  from  (5.8)  or  (5.2)  with  the  choice  of  K  used  here, 

p  »  Pyx^  =  ( l~c)ppjg  ’  (0)/(  l-*1) »  and  choosing  Er  and  E^  to 

be  negatively  correlated  results  in  a  sign-switching 

autocorrelation  for  (X  }  and  (Y  }. 

*  n'  1  nJ 

(111)  Using  Ojj  =  =  0  may  yield  a  process  with  ARNA(l.l) 

correlation  structure.  The  bivariate  system  is  now 


Px(k)  = 
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(7.6) 


Xn  *  Kn\-1  +  <>-“l2>En 

Y„  ■  ♦  <1-a22>IV 

From  (7.4)  and  (7.5)  we  have  p^l)  =  a12p  Pyy(l)  =  “22* 

where  p  =  ^(0) ,  as  before.  Further.  Equation  (7.3)  reduces 

to  p(i)  =  a^p^i-l),  t  1  2.  Hence.  Pyy(£)  =  a^.  t  i  0.  and 

PXX^^  =  (a12p)a221,  *  ^  **  Thus>  (Yn)  18  a  Process  with 
AR(1)  correlation  structure,  as  may  be  seen  directly  from 
(7.6).  In  addition,  {XR}  is  a  process  with  ARHA(l.l) 
correlation  structure,  provided  Pyx^®)  *  0.  A  particularly 
simple  way  to  achieve  this  is  to  take  =  E^.  The  process 
given  by  (7.6)  is  then  very  closely  related  to  the  EARMA(l.l) 
and  DARMA(l.l)  processes  of  Jacobs  and  Lewis  (1977,1963).  On 
the  other  hand,  we  can  choose  E^  negatively  correlated  with 
En>  e.g.,  its  antithetic,  and  may  then  Induce  negative 
correlation  in  {X^}. 

(iv)  In  order  to  generate  the  autocorrelation  function  corresponding 
to  an  AR(2)  process,  it  is  necessary  that  p(£)  satisfy  the 
difference  Equation  (7.3)  for  i  -  1  also.  Using  p(-l)  =  p(l) 
this  condition  becomes  p(l)  =  ^^/(  1— ♦g) -  Rewriting  this  with 
the  appropriate  values  for  fj,  ^2»  and  using  the  definition 
of  Pxx^^  given  by  (7.4)  yields  the  condition 

_  °22  *  all^a12a2l"<Xlla22^  (J 

P  a12(1+alla22"a12a2T5  ’ 
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where  p  *  A.(0) ,  as  before. 


In  a  model  In  which  X_  and 

n 


Yr  are  independent  the  condition  becomes 


“22  +  “l  1  ^a12°21~“l l“22^  *  °* 


(7.8) 


In  this  latter  case,  the  condition  (7.8)  may  be  satisfied  by 
choosing  either:  (a)  =  1;  or  (b)  ^  =  0  and 
alla12“21  ~  ^ase  (a)  leads  to  example  (i)  above,  i.e.,  two 
Independent  processes,  each  with  AR(1)  correlation  structure, 
one  of  which  is  degenerative,  Xr  =  X^_j .  Case  (b)  leads  to  a 
variety  of  possible  models  similar  to  those  discussed  in  the 
previous  examples,  and  one  we  have  not  noted  yet,  a  process 
with  Moving  Average  of  order  1,  MA(1).  correlation  structure. 
If  otjj  =  Ogi  =  0^2  =  0.  the  bivariate  process  is 


X„  -  KiS,Yn-l  +  <‘-“12^ 


Y  =  E' . 
n  n 


Thus,  if  E  =  E',  then  (X  }  has  the  correlation  structure 
n  n  n 

of  a  MA(1)  process  and  (Yn)  is  a  sequence  of  i.l.d. 
exponentials. 

Hence,  we  cannot  derive  an  AR(2)  structure  using  (7.8).  We  require 

dependence  between  X  and  i  and  the  satisfaction  of  condition 

n  n 

(7.7).  Since  the  relationship  between  Py^C®)  and  p^,(0)  given  hy 
(5.2)  and  the  condition  (7.7)  are  fairly  complex,  we  shall  not  attempt 
any  general  analysis.  We  note  only  that  it  is  certainly  possible  and 
give  an  example. 
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If  (E  ,E‘)  have  an  Instantaneous  correlation  of  0.5,  then 
'  n  n' 

p  =  4/17,  from  (5.2),  and  (7.7)  is  satisfied  by  =  0,  *  0.5, 

ctgj  =  0.3  and  =  0.1.  This  choice  of  K  and  the  correlation 

between  and  E^  yields  a  process  {X^}  whose  autocorrelation 

function  satisfies 


Pyx(i)  =  O.lp^f-l)  +  0.15Pxx(«-2), 


for  all  12  1*  In  passing,  we  may  note  that  (Yn)  in  this  process  has 

the  correlation  structure  of  an  ARMA(2,1)  process. 

Ve  now  consider  the  general  case,  i.e.,  the  process  with  ARHA(2,1) 

correlation  structure,  and  note  some  of  its  properties  here.  Solutions 

of  the  difference  Equation  (7.2)  depend  upon  the  behavior  of  the  roots 

o 

of  the  quadratic  equation  Z  -  ^  =  0.  It  is  easily  verified  that 

with  our  restrictions  on  {a^}  we  have:  +  ^2  ^  ®  ^  ^1  ^  ^  and 

-1  <  ^  <  1-  These  ensure  that  the  roots  of  the  quadratic  lie  within 

the  unit  circle  (See  Box  and  Jenkins,  3976,  p.  58-59,  for  a  discussion 

of  this).  However,  we  may  also  demonstrate  that  the  roots  are  real. 

2 

The  condition  for  this  is  ^  +  4^  2  0,  and 


+  4$0  =  (0tjj+a22^  + 


12  21  11  22' 


=  (aii_a22^  +  ^12^1  ^ 


Figure  7.1  shows  the  region  in  which  (^.fg)  will  lie  for  these 
processes.  It  is  worth  comparing  this  region  with  the  more  general  one 
available  for  a  stationary  process  with  an  AR(2)  component.  It  is  given 
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in  Box  and  Jenkins  (1976.  p.  59).  and  is  defined  by  +  *2  <  l’ 

f  -  f <  -l,  -1  <  f0  <  1.  In  our  case,  only  positive  is  possible, 

and  the  characteristic  quadratic  has  only  real  roots.  The  roost 
important  effect  of  this  latter  property  is  that  the  autocorrelations 
are  all  positive  if  p(l)  is  positive.  We  have  seen  in  the  examples 
above  that  this  restriction  is  not  necessary  in  the  cases  when  the 
autoregressive  component  is  first  order.  Also,  we  can  introduce 
negative  autocorrelation  at  lag  one  using  (7.4)  and  (7.5)  and  negatively 
correlated  (E  ,E').  This  may  persist  for  higher  lags  depending  on 
and  +2  as  *n  the  numerical  examples. 


Figure  7.1.  Region  of  (♦j.^  for  the  bivarlate  exponential  process. 
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8.  JOINT  DISTRIBUTION  OF  (Xn«Yn) 

We  turn  our  attention  to  the  formulation  of  the  joint  distribution 
of  (X^.Y^)  ky  investigating  the  Joint  Laplace-Stielt jes  transform. 
For  the  initial  bivariate  model  described  in  (3.1)  and  (3.2).  the 
expression  for  the  transform 


*x  Y  (s.t)  =  l{exp(-sXn-tYp)}.  (8.1) 

n  n 

is  very  complicated. 

However,  two  special  cases  are  easily  verified, 
for  all  n  in  (3.1)  and  (3.2),  we  have 
simplification 

=  [ttt]  [ttf]  -  <8-2> 


If  K<">.k£>«0 
after  considerable 


This  tallies  with  the  result  of  Section  3  that  under  these  conditions 

(Xp)  and  {Yp}  are  independent  univariate  TEAR(l)  processes,  the  first 

k  & 

with  correlation  structure  p^(k)  =  a^,  the  second  with  p99(<)  *  g^,. 
If.  in  addition.  =  Kp  and 


1  w.p.  a 

0  w.p.  1-a 


(8.3) 


then  the  distribution  of  (X  ,Y  )  turns  out  to  be  Moran's  Bivariate 

v  n  n7 

Exponential  distribution.  Using  (3.1)  and  (3.2),  we  have  from  (6.1) 


KX  Y 
n  n 


(s.t)  =  «[exp[-s{KnXn_1+(l-a)En)  -  t^Y^+O-aJE;}]]  .  (8.4) 
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Using  the  independence  of  (En)  and  {E^} ,  the  stationarity  of  {X^Y^} 
and  conditioning  on  Kn>  we  have  after  simplification 

^XyCs.t)  =  (i+s)(i+t)  -  ast  *  (8,S) 

Equation  (8.5)  is  the  Laplace-Stleltjes  transform  for  the  Moran 
Bivariate  Exponential  Distribution  as  given  in  Johnson  and  Kotz  (1970. 
p.  267). 

9.  CONCLUSIONS 

The  NEAR(2)  model  has  been  demonstrated  to  yield  a  Bivariate  First 
Order  Autoregressive  process  with  exponential  marginals  by  cross¬ 
coupling  and  auto-coupling  the  two  marginal  processes.  The  process  has 
the  same  correlation  structure  as  the  Gaussian  Bivariate  AR(1)  process. 
The  results  for  another  bivariate  process  with  exponential  marginals 
that  was  proposed  by  Kaftery  (1982)  using  the  NEAR(l)  structure  of 
Lawrance  and  Lewis  (1981)  hold  only  in  very  special  cases.  Even  then, 
as  we  have  shown,  the  correlation  structure  is  Identical  to  the  Gaussian 
AR(1)  model. 

The  possibility  of  negative  correlations  was  explored  using  the 

ideas  of  correlated  innovation  and/or  cross-correlated  attenuation. 

These  situations  occur  frequently  in  modelling  physical  phenomena.  For 

example,  the  same  shock  at  time  n  to  a  system  produces  related  effects 

in  components  given  by  X  and  Y  .  Likewise,  if  X  and  Y  are 

n  n  n  n 

flows  in  a  river  at  two  different  points  along  a  bank,  then  Er  = 
represents  a  common  phenomenon  driving  both  series.  Finally,  when  both 
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series  are  TEAR(l)  series  linked  by  a  common  attenuation  and  independent 

innovations,  the  bivariate  distribution  of  (X  ,Y  )  was  shown  to  be  the 

n  n 

familiar  Moran  distribution. 

Taken  together  these  options  Include  numerous  possibilities  to  model 
exponential  bivariate  time  series.  There  is  still  much  work  in 
parameter  estimation  before  these  models  can  be  widely  applied.  They 
certainly  lend  themselves  at  this  time  to  an  analysis  via  simulation. 

One  other  detail  which  could  extend  the  utility  of  the  model  needs 
to  be  pointed  out.  This  is  that  in  (3.1)  could  be  replaced  by 

Yn_£.  and  Xn_^  in  (3.2)  could  be  replaced  by  Xn_^,.  Th®  ®od®l  i® 
still  well  defined  and  has  obvious  physical  interpretation.  Properties 
of  this  extension  of  the  model  will  be  addressed  elsewhere. 
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